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o 

(D We investigate the properties of small-amplitude inertial waves propagating in a differ- 

Q entially rotating incompressible fluid contained in a spherical shell. For cylindrical and 

I shellular rotation proflles and in the inviscid limit, inertial waves obey a second-order 

partial differential equation of mixed type. Two kinds of inertial modes therefore exist, 

r- — I depending on whether the hyperbolic domain where characteristics propagate covers the 

Pm whole shell or not. The occurrence of these two kinds of inertial modes is examined, and we 

^ show that the range of frequencies at which inertial waves may propagate is broader than 

(~| with solid-body rotation. Using high-resolution calculations based on a spectral method, 

CL| we show that, as with solid-body rotation, singular modes with thin shear layers following 

Q short-period attractors still exist with differential rotation. They exist even in the case of 

;-H a full sphere. In the limit of vanishing viscosities, the width of the shear layers seems to 

tj^ weakly depend on the global background shear, showing a scaling in E^'^ with the Ekman 

I I number E, as in the solid-body rotation case. There also exist modes with thin detached 

layers of width scaling with E^'"^ as Ekman boundary layers. The behavior of inertial waves 

^2 with a corotation resonance within the shell is also considered. For cylindrical rotation, 

[^^ waves get dramatically absorbed at corotation. In contrast, for shellular rotation, waves 

■^ may cross a critical layer without visible absorption, and such modes can be unstable for 

cn small enough Ekman numbers. 
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1. Introduction 

The fluid flows that pervade the interior of stars and planets are most of the time 
c^ dominated by a strong background rotation. The Earth's atmosphere, oceans, or its liquid 

core are typical examples. Such flows are studied in a rotating frame, and the dominance 
of rotation arises through that of the Coriolis acceleration over all other terms in the 
momentum equation (the centrifugal acceleration is usually combined with the pressure 
gradient through the so-called reduced pressure). The Coriolis acceleration implies that 
the response of a rotating fluid to sm all ampl i tude oscillations may become singular. 

It is known indeed since the work of ICartan ( 1922 ) that the pressure perturbations of an 



inviscid incompressible fluid obey Poincare equation, which has the remarkable property 
of being spatially hyperbolic. Since boundary conditions have to be met by solutions, 
oscillations of bounded rotating fluids are a mathematically ill-posed problem. It means 
that, in most general containers, eigenmodes of inviscid incompressible fluids may simply 
not exist. Mathematically, the point spectrum of the operator is said to be empty. However, 
a few containers such as the full ellipsoid or the cylindrical layer allow the existence of 
such eigenmodes, because of the separatio n of spatial v ariables that is possible i n these 



geometries. Solutions were worked out by Kelvin (1880) for the cylinder, and by Bryan 
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( 1889 ) for the sphere and the eUipsoid. The singular nature of the solutions to Poincare 



equati on h ave been uncovered almost a century lat e r by th e works of IStewartson fc Rickard 
(1969) and Stewartson (1971 1972a). Stewartson (1972b) showed that these singularities 
sEouId turii into oscillating snear layers through, viscosity. 

New inte rest into the s i ngula r solutions to the Poincare equation was put forward with 
the work of Maas & Lam ( 1995 ) on gravity modes in two-dimensional basins . The pressure 
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fluctu ations ot internal gravity modes also satisfy Poincare equation, and^ 

(1995) showed the crucial role of characteristic trajectories and coined the concept of 
attractor t o describe the convergence of characteristics to special limit cycles. At the 
same time, HoUerbach & Kerswell (1995) and [Kerswell, (1995) investigated the internal 
shear layers spawned by the critical latitude singularity (another singularity different from 
that of the attractors) , in relation with the pre c ession flows of a spheri c al she ll 



Subsequ e nt wo rk by Rieutord & Valdettaro (1997), Rieutord et al. ( 2001h and Rieu- 
tord et al. ( 12002 ) demonstrated how attractors feature the oscillations ot slightly viscous 
rotating fluids contained in spherical shells. For very small viscosities, relevant to astro- 
physical and geophysical applications, oscillating flows are conflned within shear layers 
that close ly follow attract o rs of characteristics. In the two-dimensional limit of a spher- 



ical shell, Rieutord et al. (2002) derived analytical solutions for the structure of shear 
layers. They showed m particular that the velocity fleld of eigenmodes obeys the same 
equation as that governing the quantum states of a particle trapped in a parabolic well 
(the Schrodinger equation). 

More recently, forced oscillations in rotating flu ids have bee n examine d in the context of 
tidal interactions between two celestial bodies (Q gilvie fc Lin 2004; Wu 2005 jOgilvie fc Lin 
20071 [Goodman fc Lackner|r2009| |Ogilvie|[2009[ [Rieutord fc Valdettaro|p010p . Dissipation 
of tidally-torced modes ot oscillations m each companion alters the orbital elements and 
spin of each body. It plays a prominent role in the late evolution of close-in extrasolar 
planets. The efliciency of tidal dissipation mechanisms is usually parametrized by a single, 
dimensionless tidal quality factor Q, which is an inverse measure of the tidal torque acting 
on a gravitationally perturbed body. However, tidal dissipation is extre mely sensitive 



to the frequency and a mplitude of the external tidal potential (see, e.g., Ogilvie 2009 
Rieutord fc ValdettarOn2010) . This dependency is intimately related to the existence o 
singularities arising trom the hyperbolic nature of the Poincare equation. As in free modes 
of oscillations, these singularities consist of the critical latitude singularity and attractors 
of characteristics. 

Works examining tidal dissipation mechanisms in rotating fluids have used simple fluid 
models, where the background flow disturbed by inertial modes is in pure solid-body 
rotation. The rotation pattern inside stars and planets is much more intricate, as may 
be guessed. Large-scale flows are usually present in stars and planets. Most frequently, 
such flows are largely axisymmetric, and include differential rotation as well as a weak 
meridional circulation. In light of previous works on rigidly rotating fluids, the present 
paper is a flrst study of how departure to solid-body rotation impacts the oscillation 
properties of rotating fluids and their associated tidal dissipation. In this flrst paper, we 
focus on free inertial modes of oscillations. Tidally-forced inertial modes will be examined 
in a future study. As a flrst simpliflcation, we discard meridional circulation, which is 
very weak in general, especially for very small viscosities. However, the choice for the 
differential rotation pattern, which we assume to be axisymmetric, remains open. The 
azimuthal velocity of the (smooth) differentially rotating background flow can be expanded 
into spherical harmonics. 



M<^(r, I 



Y,^l{r)deY^ 



as for any axisymmetric azimuthal component of a vector fleld (Rieutord 1987). In this 



expression, r, 6', and (p denote the usual spherical coordinates, Y^ the spherical harmonics 
normalized on the sphere of unit radius, and wfy^{r) accounts for the radial variation in the 
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flow. If differential rotation is smooth, we expect the previous expansion to be dominated 
by small €s. In the present work, we will primarily focus on a single term, the £=l-term, 
such that the background velocity field takes the form 

Uip{r, 6) = w{r) sin 6 = rQ{r) sin 9, 

with differential rotation along the radial direction only {Q denotes the local angular 
frequency of the fluid) . The obvious asset of this functional form is that it preserves the 
relative simplicity of the coupling of harmonics imposed by the Coriolis term. Thi s type 



of flo w is known in the astrophysical literature as shellular differential rotation (Zahn 



1992). It originates from a simple description of the baroclinic flow (the thermal wind m 



geophysics), w hich arises in stably stratifled radiative zones of rotating stars (see, e.g. 



Rieutord 2006). The rotation is assumed to be uniform on spheres as a conseque nce of 



some shear-mduced turbulence, which erases latitudinal gradients of the velocity (Zahn 
1992). Low-frequency perturbations of such flows are driven by the combined effects of 
the Coriolis acceleration and buoyancy. Assuming differential rotation along the radial 
direction therefore comes to neglecting the effects of buoyancy on the fluid perturbations. 
Our model approaches this astrophysical conflguration in the limit of weak stratiflcation, 
when the Brunt- Vaisala frequency is sm all compared to the C oriolis frequency. When 



this is the case, some modes, called HI in Dintrans et al. (1999), are essentially inertial 
modes slightly perturbed by buoyancy (i.e., the Coriohs acceleration is the main restoring 
force). In addition to giving a first look at the effect of differential rotation on inertial 
oscillations, shellular rotation may also provide a (partial) view on the properties of waves 
propagating in radiative zones inside rotating stars. 

To complement the study of shellular flows, which may mimic stellar baroclinic flows, 
we also consider differential rotation associated with a barotropic fluid model. In this 
case, differential rotation only depends on the cylindrical radial coordinate s = rsin6'. 
Cylindrical differential rotation is expected in the convective zones of rapidly rotating 
stars or giant planets. It has been much investigated in the context of the Taylor-Couette 



instability, but in cylindrical geometry. H owever, recent experimental work by Kelley et al 
( [2007[ 2010J) and Rieutord et al. (2012) has raised the question of how the properties of 



inertial waves in spherical geometry are changed by the inclusion of cylindrical rotation. 
Since spherical geometry combined with cylindrical rotation results in the variables non- 
separability, which is a necessary condition for the existence of singularities like attractors, 
some properties of the perturbations can be derived by the sole analysis of the dynamics 
of characteristics. 

This paper is organized as follows. Our physical model is described in § |2j The ideal 
case of an inviscid differentially rotating fluid is flrst examined in § [3] through a linear 
analysis in the short-wavelength app roximation. Cylindrical and shellular rotation pro- 



files are considered in § 3.1| and 3.2, respectively. A detailed study of the dynamics of 



characteristics allows us to determine the range of frequencies at which inertial oscilla- 
tions may exist in differentially rotating fluids. The presence of turning surfaces, critical 
latitudes, and corotation resonances within the fluid is investigated. The full problem of 
a viscously rotating fluid is then investigated in §|4]with (linear) numerical calculations. 
Again, results with cylindrical and shellular rotation proflles are presented. Shear layers 
are shown to closely follow the attractors of characteristics obtained for an inviscid fluid 
in the sh ort- w avelength approximation. Our simulations primarily focus on axisymmetric 



modes (§ 4.2), but a f ew e xamples of non-axisymmetric modes with a corotation resonance 



are also examined (§ 4.3). Concluding remarks and future directions are drawn in § 5 



2. Inertial modes in a differentially rotating shell: physical model 

We consider a differentially rotating viscous fluid inside a spherical shell, located be- 
tween radii rjR and i? (0 < ?7 < 1). The fluid is assumed to be homogeneous, incom- 
pressible, and of constant kinematic viscosity u. We denote by u the fluid's velocity. Any 
quantity x is decomposed as x = xq + xi, where xq describes an unperturbed background 
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quantity and xi the associated disturbance (|xi| <C |a;o|). In an inertial frame, and using 
the set of spherical coordinates (r, 9, ip), the hnearized Navier-Stokes equation reads 



h ilo^ \- 2i2oez X ui + y4(ui)e<^ 



-Vpi + z/Aui, 



(2.1) 



dt dip 

with Qq = Uq/t sin 6 the fluid's angular velocity, Mq = uq ■ e;^, e^ = cos 6'er — sin 9eg is the 
unit vector along the rotation axis, and pi denotes the pressure perturbation divided by 
the reference density. Since the fluid is incompressible, the linearized continuity equation 
reads 

V-ui = 0. (2.2) 



In Eq. (2.1), A{ui) is proportional to the shear in the background flow: 

^(ui) = r sin ^ ui ■ Vl^o- (2.3) 

In the following we consider either 

• A cylindrical rotation profile, Qq = Qq{s), with s = rsin6' the radial cylindrical coor- 
dinate. This rotation profile satisfies the Proudman- Taylor theorem for an incompressible 
homogeneous fiuid in the inviscid limit, and no additional body-force needs to be assumed 
to ensure an equilibrium state for the differentially rotating background flow. 

• A shellular rotation proflle, Qq = fio{r). In this case, differential rotation is assumed 
to be maintained by some appropriate body-force f (arising for instance from Reynolds 
stresses or from a latitudinal entropy gradient), not included in our model. This body- 
force should verify V x f = —sdVlQ/dz e,^ for Taylor-Proudman theorem to be satisfled in 
this case. 



We seek solutions to Eqs. (2.1) and (2.2) in the form exp{iQ,pt + imip), where ilp denotes 
the wave angular frequency m the inertia! frame. We introduce the Ekman number E = 
z//i?^f2i-ef , where the reference angular frequency, f^j-ef , refers to ^o{R) for shellular rotation 
or t o flq js = 0) for cylindrical rotation. Dropping all subscripts and taking the curl of 
Eq. (2.1), we are to solve: 



V X (i^pU + 2f2ez x u + r sin 6{u ■ V^)e^p) = EV x Au 



V-u = 0, 



(2.4) 



where Clp = Qp + mQ is the Doppler-shifted wave frequency (that is, the wave frequency 
in the frame rotating with the fluid). 

We complete Eqs. (2.4) with stress-free boundary conditions at both boundarie s: u-er = 
g]er = 0, where \a] denotes the viscous stress tensor. As discussed in iRieutord 



and e_c^ x 



& Valdettaro (2010), the general properties of oscillation modes are essentially unchanged 
by the choice ot boundary conditions. We thus opt for stress-free boundaries, which are 
numerically less demanding. 



Inviscid problem: paths of characteristics, existence of turning surfaces and 
corotation resonances 



Before examining the full solutions to Eqs. (2.4) by means of numerical simulations, 
we flrst study the inviscid limit and especially the dynamics of characteristics, which is 
known to feature the eigenfunctions of the viscous problem. For this purpose, we shall use 
cylindrical coordinates {s,ip,z), and eliminate the velocity perturbations in favor of the 
pressure perturbations ( p). T he partial differential equation (PDE) satisfle d by p is given 
in the appendix by Eq. (A 6). Retaining only the second-order terms, Eq. (A6) becomes: 



d'^p Az d'^p 
ds'^ d"^ dsdz 



+ 1 



d'^p 
dz'^ 



0, 



(3.1) 
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where 



A.. 



2Q d 
s ds 



[s^Q) and A^ 



2Q d 

s dz 



{s^n) . 



(3.2) 



Eq. (3.1) matches equation (A 16) of [Ogilvie (2005) for an homogeneous incompressible 
fluid. The equation governing the paths ot characteristics in a meridional plane then reads 



dz 
ds 



2^2 



A,±i 



1/2 



with e(s, z) = Al + ml{A, - nl 



(3.3) 



For sohd-body rotation, A^ = 0, As = 4^2^^ and Eq. (3.1) reduces to the well-known 
Poincare equation for inertial waves. In this case, characteristics form straight lines in a 
meridional plane, and the fact that ^ is constant makes the entire shell either hyperbolic 
or elliptic for inertial waves. With differential rotation, however, the dependence of Ag, 
Az and f2 with s and z implies that characteristics no longer form straight lines, as will 
be illustrated in the following. Also, ^ may vanish and change sign inside the shell. This 
makes possible the co-existence of hyperbolic and elliptic domains for pure inertial waves 
in a differentially rotating spherical shell, separated by turning surfaces. It is reminis- 
cent of the existence of turning surfaces for gravito-inertial waves and magneto-inertial 
waves in rigidly rotating spherical shells, for which the p erturbed pressure also satisfies 
a second-order PDE of mixed type in the in viscid limit (iFr iedlander fc Siegmann 1982 



Friedlander 



1982, 1987 Dintrans et al. 1999). Equation (3.3) also shows that, in contrast 



to solid- body rotation, paths ot characteristics with differential rotation depend on the 
azimuthal wavenumber m through fip. In the rest of this section, we study in more de- 
tail the dynamics of characteristics based on Eq. (3.3), considering a particular profile of 
cylindrical rotation in § 3.1 , and of shellular rotation in § 3.2 



3.1. Cylindrical rotation 

3.1.1. Paths of characteristics 

To restrict the differential rotation function space to a one-parameter space, we consider 
the following background rotation profile: 



n{s)/n 



ref 



(3.4) 

with r^j-ef = ^{s = 0) the shell's angular freque ncy at the rotation axis, and where \e\ can 
be seen as a Rossby number for the fiuid. Eq. (3.1) then takes the form: 



d'^p 
ds"^ 



+ 1 



d'^p 
dz^ 



0, 



(3.5) 



where 



fip(s) 



is the Doppler-shifted frequency, and 

As{s) =4fi2(s) X 



i7p + mi7(s) 
, <slKf 



1 + e^sjKf 



(3.6) 



(3.7) 



We will take e > -1/2 so that Ag > Vs G [0, R] (Rayleigh's stability criterion). We will 
denote the quantity As by k^ (it is known as t he s quare of the radial epicyclic frequency 
in astrophysical discs). We point out that Eq. (3.5) is symmetric by z — >■ —z, so that only 
positive values of z will be considered. From Eq. ( 



dz 
ds 



±e 



1/2 



3.3), paths of characteristics satisfy 
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with 

^(s) = -^ - 1. (3.9) 

3.1.2. Existence of turning surfaces: D and DT inertial modes 

The relation C,{s) = defines turning surfaces that separate the hyperbohc and elhptic 
domains. They correspond to cyhnders in a meridional slice of the shell. Their location 
satisfies 

Qlis) = kI{s), (3.10) 

with waves propagating when ^^(s) ^ ^p{s). Here and in the following, we will assume 
that the hyperbolic domains of the spherical shell, where ^ is positive, host inertial modes 
of oscillations. The dependence of .^ on s leads to two kinds of inertial modes with cylin- 
drical rotation: 

(i) modes that exhibit at least one turning surface within the shell, which we call DT 
modes for future reference (D for differential rotation, and T for turning surface), 

(ii) generalized inertial modes with no turning surfaces within the shell, which we name 
D modes. 

The occurrence of axisymmetric {m = 0) and of some non-axisymmetric D and DT 
modes is illustrated in Fig. [T] In all panels, the mode's eigenfrequency in the inertial 
frame {fip) is shown in x-axis, and the differential rotation parameter e is in y-axis. 



Contours are obtained by numerically calculating the right-hand side of Eq. (3.9). DT 
modes, which are represented by white areas, have ^ < at least once within the shell. D 
modes satisfy ^ ^ everywhere in the shell and are shown by the blue area. No inertial 
modes can propagate if ^ < in the entire shell, such case being depicted in red. Note that 
axisymmetric inertial waves may propagate for eigenfrequencies Qp > 2r2i.ef if e > (see 
top-left panel in Fig. [1]). Note also that the results shown in Fig. [llare independent of the 
shell's aspect ratio i], as is actually the dynamics of characteristics with our cylindrical 
rotation profile. 

The transition curves that separate D and DT modes satisfy ^(s = 0) = and 
^{s = R) = 0. The solution to ^(s) = can be cast as fip/fij-ef = —m[l + e{s/ Kf'] ± 
2a/[1 + 2e{s/RY][l + e{s/RY]. Applied to the following cases: 

• ^(s = 0) = 0, that is when the turning surface matches the location of the rotation 
axis, we have 

fipMef = -"^±2, (3.11) 

such case being depicted by vertical dotted lines in Fig. [ij 

• ^(s = i?) = 0, when the turning surface passes by the shell's surface, yields 



fipMef = -m{l + e)± 2^(1 + 2£)(1 + e), (3.12) 

which is overplotted by solid curves in Fig. [l] 

We point out that there may exist two turning surfaces in a meridional plane. This occurs 
for instance for the DT modes with m = A located in the (small) region of parameter 
space below the solid curve near Vtp = —'^^^rei and e = —0.5. 

Non-axisymmetric modes may feature a corotation resonance within the shell when 
flp{s) = 0. The location of such critical cylinders is given by 



.„i. = fl,/.-'x(^-l--A-). (3,13) 

For m > 0, modes with £ ^ have a critical layer if —m ^ Qp/Qj-^^ ^ —m{l + e), and 
those with e ^ if —m{l + e) ^ Qp/Qj-g.^ ^ —m. Modes that exhibit a critical cylinder 
within the shell are located between the two orange solid curves in Fig. [l} 
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Figure 1. Illustration of the two kinds of inertial modes propagating in a differentially rotating fluid, with 
cylindrical rotation profile Q,{s)/Q,rct = 1 + £{s/ Rf' (f2rcf denotes the shell's angular frequency at the rotation 
axis). Results are shown for azimuthal wavenumbers m equal to 0, 1, 2 and 4 from top-left to bottom-right. DT 
modes, represented by white areas, exhibit turning surfaces inside the shell. D modes correspond to the blue 
area and do not feature turning surfaces inside the shell. No inertial modes may propagate in the red areas. The 
solutions to 5(s = 0) = are shown by vertical dotted lines, and the solutions to ^(s = R) = Q hy solid black 
curves. Non-axisymmetric modes that feature a critical cylinder inside the shell (corotation resonance) are located 
between the two orange solid curves. 



The minimum and maximum frequencies at which inertial waves may propagate in a 
shell with cylindrical rotation is summarized in Fig. |2] for the four azimuthal wavenumbers 
considered in Fig. fl} These frequencies are obtained oy solving ^{s) = 0. They corre spond 
to th e min imum and maximum frequencies amongst the four values given by Eqs. (3.11) 
and (3.12), except for m = 4 and e approaching —0.5 (in which case, as already pomted 
out before, there are two turning surfaces inside the shell). Fig. |2] highlights that the range 
of eigenfrequencies at which inertial modes of oscillations may exist in a spherical shell 
is broader with cylindrical rotation than with solid rotation, particularly if the rotation 
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Figure 2. Range of eigenfrequencies (fip, x-axis) for inertial modes of oscillations in a rotating shell with cylindrical 
rotation profile Q,{s)/Q.y.ct = 1 + e{s/ RY . Eigenfrequencies are in units of firef, the shell's angular frequency at 
the rotation axis. Results are for m = 0, 1, 2 and 4. For a given m, the minimum and maximum frequencies at 
which inertial modes exist correspond to the left and right curves, respectively. 

profile is an increasing function of the cylindrical radial coordinate. This is expected, 
because of the larger frequencies at which the fluid may rotate in this case. 

3.1.3. Dispersion relation, phase and group velocities 

We examine in this paragraph the propagation properties of local waves in a meridional 
plane. In the short-wavelength approximation, the general expressions for the wav e dis - 
persi on rela tion, the phase and group velocities are given i n the Appendix by Eqs. (A 7), 
(A 9) and (A 10). For our cylindrical rotation profile, Eq. (A 7) reduces to 



n^is) 



ki 



12' 



(3.14) 



where ||k|| = y fcg + ^z- I^ the frame rotating with the fluid, the phase velocity, Vp 

fCrK 



fipk/llklP, satisfies 



Vn = ±K,Js) 



\W 



while the group velocity, Vg = dClp/dk, is given by 



k 
Vg = it/ts(sj--—j -g (— /c^jeg + Kge^.)- 



(3.15) 



(3.16) 



From Eqs. (3.9) and (3.14), fcg = at a turning surface (where ^ vanishes). At this location, 
the radial component ot the phase velocity therefore vanishes, as well as the two other 
components of the group velocity (that is, rays are perpendicular to the turning surface). 

A corotation resonance occurs inside the shell where ^p{s) = 0. This formally happens 
if l^sl — )■ oo, in which case the modulus of the phase velocity and the group velocity both 
tend to zero, or if A;^ — !■ 0, in which case the phase velocity cancels out while the radial 
component of the group velocity cancels out. In both cases, characteristics become more 
and more parallel to the rotation axis as they approach the corotation resonance (their 
slope \dz/ds\ increases), and they do not cross corotation (as the radial component of the 
group velocity tends to zero at this location). The fact that the phase velocity progressively 
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decreases to zero upon approaching corotation suggests that non-hn ear effects should 
become important near this location, as pointed out for instance by Barker fc Ogilvie 
(2010) for internal gravity wa ves in a background shear flow, which satisfy a dispersion 
relation similar to that in Eq. (3.14) with the epicyclic frequency squared k^ replaced by 



the Brunt- Vaisala frequency squared. 



3.2. Shellular rotation 

We now come to the propagation properties of inertial modes of oscillation in fluids with 
shellular rotation (the background rotation profile being a prescribed function of the 
spherical radius). As in the case of cylindrical rotation, we show in this section that two 
kinds of inertia! modes can be distinguished, based on the presence or not of turning 
surfaces within the shell. We then discuss the range of eigenfrequencies at which inertial 
modes may exist, the existence of critical latitudes and corotation resonances. 



3.2.1. D and DT inertial modes 

We adopt the same strategy as in § 3.1 for cylindrical rotation and restrict the differential 
rotation function space to a one-parameter space. We take the following shellular rotation 
profile: 

fi(r) = firef X (r/i?)^ (3.17) 



where r = vs^~+^, fij-ef = fi(r = R) is the angular frequency at the shell's surface, and 
a > —2 to avoid centrifugal instabilities (Rayleigh's stability criterion, Ag > 0, is satisfied). 
The quantity |1— r^'^l, which measures the relative difference of angular frequencies between 
the inner core and the shell's surface, can be seen as a Rossby number for the fluid. 



Eq. (3.1) being symmetric by ^ — )■ —z for shellular rotation p rofil e, only positive values of 
z will be considered (as with cylindrical rotation). From Eq. (3.3), paths of characteristics 
satisfy 



where 



as, 



A' 



dz 
ds 



s2^2 



sz 



Mr)f,±e'\ 



w^^^^'-^w 






A{r) 



i?2 Qyr) 



a 



^lir) 



(3.18) 



(3.19) 



(3.20) 



Vtp + mVL{r) is the Doppler-shifted wave frequency. Whenever paths of char- 



and VLp{r) 

acteristics are calculated by solving Eq. (3.18) (or Eq. (3.8) for cylindrical rotation), 
reflections are imposed at the shell's inner core (r = rjR) and surface (r = i?), at the 
rotation (s = 0) and equatorial {z = 0) axes, and at the location of turning surfaces 
(where ^ = 0). 

The possibility that ^ vanishes within the shell leads again to two kinds of inertial modes 
with shellular rotation: DT modes, with at least one turning surface inside the shell, and 
D modes, with no turning surfaces inside the shell. The occurrence of axisymmetric and 
of some non-axisymmetr ic D and DT modes is depicted in Fig. [3j Contours are obtained 
upon calculation of Eq. (3.19) in a spherical shell with t] = 0.35. As in Fig. [11 D and DT 
modes are depicted by blue and white areas, respectively, and no inertial modes exist in 
the red areas. 

The relation ^ = defines turning surfaces with colatitude 6{r) given by 



sm 



1 + 



2ni(r) 



a 



fi2(, 



±li'l\r) 



(3.21) 
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Figure 3. Occurrence of D and DT inertial modes with the shellular rotation profile i}{r)/Q,,.cf = {rjR)" (f^ref 
denotes the angular frequency at the shell's surface). Results are obtained for r) — 0.35, and m = 0, 1, 2 and 4 
from top- left to bottom-right. DT modes (white areas) have at least one turning surface inside the shell. D modes 
(blue area) do not have turning surfaces. No i nertia l m odes exist within the red areas. The dotted, dashed, solid 
and dash-dotted curves are solutions to Eqs. (3.231 to 13.261, respectively. In the bottom-right panel, grey areas 
depict DT modes located between two turning surfaces that englobe the shell's inner core and surface. For m 7^ 0, 
modes that feature a critical layer inside the shell (where Qp{r) — 0) are located between the two orange solid 
curves. 



with 



fi{r) 



4(4 + a)nl{r) 



(3.22) 



Another way to find the occurrence of DT modes is to require that, for given values of 
m, Up and a, the inequahty ^ sin 9{r) ^ 1 is satisfied at least once within the shell. 
We also point out that the two propositions V{s, z}, ^(s, z) > 0, and Vr, |^p(r')| < 2f2(r), 
seem to be equivalent numerically (as we have checked). In other words, if the whole shell 
satisfies the condition for hyperbolicity, it also meets the standard criterion for propagation 
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Figure 4. Minimum and maximum eigenfrequencies (x-axis) at which inertial modes of oscillation may exist with 
shellular rotation profile Q,{r)/^^Ei = {I'/RVy for m = 0, f, 2 and 4 {rj = 0.35). Eigenfrequencies are in units 
of the angular frequency at the shell's surface. The inset plot focuses on previous range of eigenfrequencies for 
a & [-0.4,0.4]. 

of inertial waves everywhere locally. Thus, regions with no mode propagation, and those 
with D modes, can be distinguished by requiring that D modes verify |fip(r)| < 2Q{r) in 
the whole shell. 

The transition curves that separate D and DT modes then satisfy: 

• sin 9{R) = 0, when the turning surface hits the shell's surface at colatitude 9 = 0. 
It reads \np{R)\ =2Q{R), or 

^p/^ref = -m ± 2, (3.23) 

expression that is depicted by vertical dotted lines in Fig. |3| 

• sin 9{R) = 1, which yields 



m 



lip/ i ij.gf - 

which is overp lotted by a dashed curve. 
• sin 9{r]R) = 0, which can be recast as 

Qp/Qj-Qi = (2 — m)ri'^ 

^p/^ref = -{2 + m)r]'^ 
and are shown by solid lines. 



±^2(2 + a), 

if fip/(2 - m) > 
if Qp/{2 + m) <0, 



(3.24) 



(3.25) 



sin 9{r]R) = 1, when the turning surface hits the shell's inner core at colatitude 



9 = 7r/2, that is when 



Qp/Qrei = V"" (-'m ± 2y^l + (T/2') 



(3.26) 



expression that is shown by dash-dotted curves. 

DT modes with \m\ ^ 2 are located between the turning surface and the shell's surface 
for (T > (since the background velocity increases with radius), and between the shell's 
inner core and the turning surface for a < 0. A more complex situation arises when 
\m\ > 2, where the equation \ Qp{ri R)\ = 2Q{riR) may admit two solutions for a given 
fip, as can be seen from Eq. (3.25). In the bottom-right panel of Fig. [3] (?7i = 4), DT 
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Figure 5. Location of turning surfaces for m = modes with varying the power-law index of the shellular 
rotation profile, for ftp = l.SJlrcf (left panel), Q,p = 2r2rcf (middle panel), and fip — 2.5ilrci (right panel). 



modes located between the two solid lines propagate between the shell's inner core and 
the turning surface, since as D modes they satisfy \Clp{T]R)\ < 2Q{'r]R). In contrast, the 
DT modes above the upper solid curve, and those below the lower solid curve, propagate 
between the turning surface and the shell's surface. Interestingly, modes with \m\ > 2 
may exhibit two turning surfaces that are located between the shell's inner core and 
the shell's surface. This brings the possibility of waves naturally bouncing between two 
turning surfaces without interacting with the shell's boundaries. Such modes are depicted 
by grey areas in the bottom-right panel in Fig. [3j Note, however, that this particular set 
of non-axisymmetric modes will necessarily feature a critical layer (corotation resonance) 

between both turning surfaces, where Clp vanishes. The location of such critical layers is 
given by r = i?(— fip/TTifij-ef)^''^. For m > 0, modes with cr ^ have a critical layer if 
—mrj" ^ fip/fij-ef ^ — m, and those with a ^ if —m ^ fip/fij-ef ^ —mr]'^. Modes that 
exhibit a critical layer within the shell are located between the tw o orange solid curves in 
Fig. [3j They will be further discussed in Sections |3.2.3| and |4.3.2 

The range of eigenfrequencies Qp for which inertial modes of oscillation may exist in 
a spherical shell with t] = 0.35 is displayed in Fig. |4]for m = 0, 1, 2, and 4. Compared 
to solid-body rotation, inertial modes with shellular rotation may exist at much larger 
frequencies, especially if the rotation profile is a decreas ing fu nction of radius. 

Finally, the location of turning surfaces, given by Eq. (3.21), is depicted in a meridional 
slice in Fig. [5]for m = modes with a ranging from -2 to 2, and increasing by 0.2. Results 
are shown for Qp/Qy-^i = 1.5, 2, and 2.5 from left to right. For Qp = 2r2i.ef, turning surfaces 
exist for all values of a but for a = (corresponding to solid-body rotation; note however 
that turning surfaces exist in the limit |(t| — )■ 0). In all three panels, DT modes with a < 
exist between the shell's inner core and the turning surface (they may therefore exist 
without reflexions on the shell's surface). Similarly, DT modes with a > exist between 
the turning surface and the shell's surface, and may theref ore ex ist without reflexions on 
the shell's inner core. This case will be illustrated in Section 4.2.2 where we will show that 
DT modes may feature shear layers following attractors of characteristics in the absence 
of an inner core (the necessary reflexion being provided by the shell's surface and the 
turning surface). 

3.2.2. Critical latitudes 

Rays tangent to the shell's inner core defi ne a c riti cal la titude 6[. They satisfy dz/ds = 

—s/z with rjR = vs^+z^. Combining Eqs. (3.18) to (3.20), rays intersect the shell's inner 
core at s = Sr, with 



'c, 



r]R 1 



\ 1/2 

AQ^{r]R) j 



(3.27) 
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A critical latitude 9i at the shell's inner core therefore exists for modes with \VLp{riR)\ ^ 
2Q{r]R), and verifies 



sin 6*1 



(3.28) 



2n{r]R) ' 

which is a straightforward generalization of the expression for solid-body rotation. For 
fixed values of rj, m and Vtp, and for cr > 0, sin 6*1 —)■ 1 as cr increases. For o" < 0, sin ^j — )■ 
as |o"| increases. 

Similarly, it is easy to show that a critical latitude 9^ exists at the shell's surface for 
those modes that satisfy |i7p(i?)| ^ 2i7(i?), with expression given by 



sin^r 



2n{R)' 



(3.29) 



independently of a. The outer critical latitude has therefore the same value as in the 
solid-body rotation case. 

Critical latitudes also exist with cylindrical rotation. For o ur cylindrical rotation profile, 
an inner critical latitude exists if ^(s), given by Eq. (3.9), is positive for s G [0,riR]. 
Similarly, an outer critical latitude exists if ^ ^ Vs G [^,R] (that is, for D modes only). 
There is, however, no simple expressions for the inner and outer critical latitudes with our 
cylindrical rotation profile. 

3.2.3. Dispersion relation, phase and group velocities 



In the short-wavelength approximation, the wave dispersion relation given by Eq. (A 7) 
becomes 



^lir) 



AQ'^(r) 



kt 



X 




1- 



ZrCg 
SKz 



(3.30) 



for our shellular rotation profile, with ||k|| = y A;^ -|- A;|. We see that the bracket term 
in the above equation can be negative: inertial waves may become unstable for shellular 
rotation. As discussed in the appendix, this is a form of t he Goldreich-Schubert-Fricke 



instability i n the limit of an inviscid incompressible fluid (Goldreich & Schubert 1967 



Fricke||1968 
whenever c) 



In that limit, Eq. (A 7) indicates that the instabihty can m prmciple occur 



l./dz 7^ (or, equivalent ly, when Az 7^ 0). Assuming that inertial waves remain 



stable against our shellular rotation profile, we introduce 



B 




zkg 

SKz 



1/2 



(3.31) 



and Eq. ( 3.30[ ) reduces to fip(r) = ±2-Bn(r)|fc^|/||k||. Note that i3 = 1 for solid-body 
rotation. In the rotating frame, the phase velocity reads 



and the group velocity is given by 



±2Bn{r] 



-KzGs 1 fCs^z 



kzii- 



B 



(3.32) 



a sz ||k|p 
4 r^ kakz 



B 



-1 



(3.33) 



We now examine the case where a critical layer forms inside the shell, where ^p{r) 



0. 



Inspection of Eq. (3.30 ) shows that there are formally three cases where f2p(r) may vanish: 



fc, — )• 00 at finite k 



In this case, Vp 



and Vjr = at corotation, which is similar 



to the behavior obtained with cylindrical rotation (inertial waves do not cross corotation 
in such case). 
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• kz -^ at finite ks, in which case Vp = 0, Vg-eg = 0, and |vg-ez| — )■ oo. Surprisingly, 
this case indicates that inertial waves may propagate across a critical layer (with locally 
vertical paths of characteristics). 

• i3 — )■ 0, which implies that Vp = 0, |vg ■ eg] — )■ cxd and |vg ■ ez| — )■ oo. This case also 
formally shows that inertial waves may cross a critical layer. The condition B = can be 
recast as 



dz 
ds 



Ka 



s 

crit fcz 



= --fl + -4y (3-34) 



which shows that, if S = at corotation, then paths of characteristics have a finite slope 
upon crossing the critical layer. 

The propagation proper ties of inertial waves in viscous fluids with a critical layer will be 



examined in Section |4.3.2| with the results of numerical calculations. We stress again that 
the fact the phase speed vanishes at corotation indicates that non-linear effects should 
occur near this location, whic h are n ot captured by the present linear analysis, nor by the 
numerical solutions of Section 14.3.21 

4. Viscous problem: shear layers, behaviour at corotation resonances, comparison 
to inviscid analysis 

We have described in Section[3]the propagation properties of inertial waves in a spherical 
shell with (specific) cylindrical and shellular rotation profiles. We have shown that in both 
cases, modes of oscillations satisfy in the inviscid limit a second-order PDE of mixed type. 
Two families of eigenmodes can therefore be distinguished: D modes, which may propagate 
throughout the entire shell (the hyperbolic domain covers the whole shell), and DT modes, 
which feature turning surfaces inside the shell (the hyperbolic domain covers part of the 
shell). Compared to solid-body rotation, the presence of turning surfaces with differential 
rotation tends to broaden the range of frequencies at which inertial waves may propagate 
in a spherical shell. 

An interesting result from our analysis in Section |3] is that paths of characteristics 
depend on the azimuthal wavenumber m (through the space-varying Doppler-shifted fre- 
quency). In contrast to solid-body rotation, the structure of non-axisymmetric modes 
cannot be simply inferred from that of axisymmetric modes. However, as far as their dy- 
namical properties are concerned, m 7^ and m = modes only differ by the possible 
presence of corotation resonances for m, 7^ 0. We will therefore examine the general p rop- 
erties of viscous modes of oscillations for m = modes only. This is done in Section |4.2 
where we take the same rotation profiles as in Section |3} We show that the structure o: 
shear layers in the viscous problem closely follows the paths of characteristics derived in 
the inviscid problem. We particularly focus on singular modes featuring a shear layer that 
converges towards a short-period attractor. Some of these singular modes form through 
reflexions at turning surfaces and may therefore exist in the absence of an inner core. A 
qualitative analysis of the occurrence of modes with quasi-periodic orbits of characteris- 
tics is also done at Section 4.2.3 We finally examine in Section 4.3| a few ?7i 7^ viscous 



eigenmodes with corotation resonances. 

4.1. Numerical method 

Our numerical simulations solve the linearized differential system (2.4) and related stress- 
free boundary condit ions, using a sp ectral approach. Fields are decomposed into spherical 



harmonics (see, e.g., Rieutord 1987), with the perturbed velocity expanded as 



00 e 
u{r,e,^) = J2J2 <ir)^T + vL{r)SY^ + wUr)TY^, (4.1) 

with R™ = Yl^{e,if)ey,, S™ = VY^, and T™ = V x R™, and where Y^^ denote the 
usual spherical harmonics normalized on the sphere of unit radius. In the radial direction. 
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equations are discretized on Gauss-Lobatto collocation nodes associated with Chebyshev 
polynomials. Equations are truncated at order L for the spherical harmonics basis, and at 
order A^^ for the Chebyshev basis. The incomplete Arnoldi- Chebyshev algorithm is used to 
compute pairs of complex eigenvalues (flp) and eigenvectors, given an initial guess value 
for flp (for details on this numerical method, see Valdettaro et a/.||2007 ). Values of Nr and 
L are mode-dependent, and will be specified below. Unless otherwise stated, all modes 
are calculated assuming symmetry with respect to the equatorial plane. 



4.1.1. Shellular rotation 

For the shellular rotation profile that we consider, Vt{r)/VL^(,f 
Eqs. (2.4) onto R™ and T™ gives the following two equations: 



[r/RY, projecting 
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dr 



dr 



dQ 



+r— [-{i+l)A{i + l 



EAfAf(ruL) + 
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1+2. .1+1^ 
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2im^{ut, + V, 
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m)w, 



l-i 



+ B{i + l,m)w^+^] 



where 

A{i,m) 



m 



1/2 



£2 \ 4£2 



B{i, m) = t{t - 1)A(£, m), A^ 



(4.2) 



1 d'^ _ £{£ + 1) 



In Eqs. (4.2), since the fluid is incompressible, 

1 d{r^ui,) 



.1 



+ l)r dr 



(4.3) 



For solid-body rot ation, Vt{r) = i^^ef and Eqs. (4.2) reduce to equations (2.2) of Rieutord & 
Valdettaro (|l997 ) (note that these authors defined the Ekman number as E = i7/2K^l^, 
while our cletimtion does not include the factor 2). After projection on the spherical 
harmonics, stress-free boundary conditions read 



u^ 



dr"^ 








(4.4) 



for the radial functions taken at r = rjR or r = R. Eqs. (4.2) to (4.4) can be recast as 
an eigenvalue problem featuring a tridiagonal block ma trix, hke for solid-body rotation. 
More details into the numerical scheme can be found in Rieutord & Valdettaro (|1997j). 
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4.1.2. Cylindrical rotation 

For our cylindrical rotation profile, D.{s)/Vlj-Qf = 1 +e{s/Kf', projecting Eqs. (2.4) onto 
R™ and T^ gives, after some tedious algebra, two long differential equations, which we 
do not write for the sake of legibility. The second-order differential equation satisfied by 
w^ can be cast as a linear combination of m^^, du^^/dr, u^^, dufy^^/dr and w^^ (the 
term in front of wf^'^ vanishes for m = 0). The fourth-order differential equation for uf^ 
exhibits terms in w^"^, dwf^^/dr, w^\ dwf^^/dr, m^^, du^'^/dr and d'^u^^/dr'^ (the 

terms in front of "U^^, du^'^ /dr and d'^u^'^ /dr"^ vanish for m = 0). Along with stress- 
free boundary conditions, these equations take the form of an eigenvalue problem where 
the matrix to invert has blocks composed of (up to) seven bands. The coupling between 
harmonics of ranks i, i ±1, i ±2 and i ±3 arises from our particular cylindrical profile. 

4.2. Axisymmetric eigenmodes: a few illustrative cases 
We present in this section the results of simulations for axisymmetric eigenmodes with 
moderate viscosity, the default Ekman number being E = 10^^. Unless otherwise stated, 
our simulations are carried out with a shell's aspect ratio rj = 0.35. The propagation 
properties of shear layers for weakly-damped eigenmodes is compared with the analytic 
description based on paths of characteristics, detailed in Secti on pi Cylindrical and shel- 



lular rotation profiles are examined in Sections 4.2.1 and 4.2.2, respectively. A qualitative 
analysi s of th e occurrence of modes with quasi-periodic orbits of characteristics follows in 
Section KTSi 

4.2.1. Cylindrical rotation 

We describe in this paragraph three representative m = eigenmodes obtained with our 
cylindrical rotation profile: a D mode (the hyperbolic domain covers the whole shell), and 
two DT modes (the hyperbolic domain covers a fraction of the shell). All three calculations 
have a spectral resolution N^ = 450 x L = 800. 

Fig. |6] displays the results for a D mode with eigenfrequency Qp ~ O.Olfirgf, obtained for 
e = —0.3 (recall that for cylindrical rotation, fij-ef denotes the shell's angular frequency at 
the rotation axis). The top-left panel depicts contours of the mode's kinetic energy in a 
meridional quarter-plan^ R esult s are shown only for positive values of z, for reasons of 



symmetry exposed in Section 3.1[ This D mode essentially features a shear layer following 



a short-period wave attractor, which is overplotted by a thick white curve. The mode's 
kinetic energy is maximum ne ar the rotation axis, like in the case of solid-body rotation 



(Rieutord & Valdettaro 1997). Note also the presence of another shear layer tangent to 
the shell's inner core. As already highlighted in Section [3TT| rays in a meridional plane are 



curved by differential rotation. This is also illustrated in the t op-r ight panel of Fig.pl which 



displays the paths of characteristics obtained by solving Eq. (3.8 ), given a particular initial 
location (here at s = and z = 0.75R). The patterns of the viscous shear layer (left panel) 
and of the paths of characteristics under the short-wavelength approximation (right panel), 
are in excellent agreement. Paths of characteristics focus towards the aforementioned wave 
attractor, depicted by a red curve in this panel. Upon varying initial conditions, we have 
checked that this is the only attractor in the shell for this particular eigenfrequency. The 
bottom panel of Fig. |6]shows the good convergence in spectral resolution for this mode. 
We present in FigTrT] the results for an m = DT mode with Vlp ^ 1.33^1^^^ and 
e = —0.4. The left panel shows a meridional cut of the kinetic energy for this mode, 
which exhibits again a shear layer following an attractor of characteristics. The mode's 
pattern features reflexions on the inner and outer edges of the shell, as well as on a turning 
surface located at s ~ 0.75i?. The radius at which the shear layer bounces off is in good 

f In all our contour plots, the simulation parameters are summarized below the x-axis (spectral resolution 
Nr X L, azimuthal wavenumber m, Ekman number E, shell's aspect ratio i], differential rotation rate e or a for 
cylindrical or shellular rotation profiles, respectively). Similarly, the complex eigenfrequency is indicated in the 
bottom-left corner of the shell: the absolute value of its real part is denoted by |aj| (in units of firet), and its 
imaginary part by r (in units of fl^^]; negative values of r correspond to damped modes). 
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Figure 6. Top left: meridional cut of the kinetic energy for a D mode with eigenfrequency fip ~ 0.91f2rcf, 
obtained with E — 10~* and cylindrical rotation profile Q{s)/Q,,-c{ — 1 — 0.3(s/i?)^. An attractor of characteristics 
is overplotted by a solid white curve. Top right: a few paths of characteristics are shown to converge towards the 
aforementioned attractor, which is displayed here by a red curve. Bottom: spectral content of the radial (u) and 
orthoradial (w) components of the velocity field for this mode. Chebyshev and spherical harmonic coefficients are 
shown in the left and right panels, respectively. 



agreement with the (numerical) solution to ^(s) = 
of characteristics, obtained again by solving Eq. (3.8|) 



shown by a dashed line. A few paths 
are shown by white curves in the 
right panel of Fig. [tI which rapidly converge toward the above attractor (in red in this 
panel). We point one that there exists another attractor for this eigenfrequency and this 
differential rotation rate, which is overplotted by a blue curve. 

An interesting DT mode is finally illustrated in Fig. [8] for Qp ^ LlOfij-gf and e = —0.4, 
where a shear layer follows an attractor for an arbitrarily small inner core. In this particular 
calculation, rj = 0.05. This example highlights that, in contrast to solid-body rotation, 
attractors may be formed by reflexions on one of the shell's edges (the surface in this 
case) and a turning surface. T his property is also shared by gravito -inertial modes in 



stably-stratified rotating fluids ( Friedlander||1982 Dintrans et al. 1999). The mode shown 



in Fig. [8] indicates that, with cylindrical rotation, wave attractors may exist independently 
of a putative inner core inside the shell. We show below that the same is true with shellular 
rotation. 

4.2.2. Shellular rotation 

We describe in this paragraph a few axisymmetric {m = 0) eigenmodes obtained with 
our shellular rotation profile, f2(r)/fii-ef = {r /RY (here, fij-ef denotes the angular frequency 
at the shell's surface) . An overview of the impact of shellular rotation is illustrated in Fig. |9] 
for a singular mode with eigenfrequency Vtp f» Lllfij-ef for solid-body rotation. Contours of 
the kinetic energy in a meridional quarter-plane are shown for a equal to (top-left panel) , 
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Figure 7. Left: meridional slice of the kinetic energy for a DT mode with Qp « 1.33f2rcf, E — 10~ , and e = —0.4. 
The shear layer fo llow s a wave attractor (solid curve). The dashed line shows the location of the turning surface 
inferred from Eq. (3.9 1. Right: a few paths of characteristics are depicted by white curves, which converge toward 
the aforementioned wave attractor (red curve). A second wave attractor is shown in blue. 
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Figure 8. Meridional slice of the kinetic energy for a DT mode with Up ^ l.lOQrct, E = 10~ , 77 = 0.05 and 
e = —0.4. The turning surface is shown by a dashed line. For an arbitrarily small inner core, the shear layer follows 
an attractor of characteristics, shown by the solid curve. 



—0.5 (top-right), and —1.3 (bottom- left). Calculations were carried out as a sequence, the 
eigenfrequency of the least-damped mode obtained with differential rotation rate a being 
used as the guess frequency for the next run with differential rotation rate a — da. In this 
series of runs, da = 0.01, E = 10~^, and the spectral resolution is Nr = 250 x L = 600. 
The contour panels in Fig. |9] show how the shape of the shear layer, and thus that of the 
wave attractor it focuses on (which curve), is progressively deformed by shellular rotation. 
Note the presence of a turning surface for a = —1.3 (dashed curve). Not surprisingly, the 
mode's eigenfrequency increases with increasing the steepness a of the rotation profile. 
We also observe that the mode's damping rate decreases in absolute value for steeper 
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Figure 9. Meridional cuts of kinetic energy for a given eigenmode with E = 10~^ and a = (top-left), —0.5 
(top-right), and —1.3 (bottom-left). The shear layer for this mode follows an attractor overplotted by a solid curve 
in each of the se panels. The dashed curve in the bottom-left panel displays the location of a turning surface, given 
by Eq. (3.21 1. The bottom-right panel shows that the damping rate of previous mode with a = —I scales as E^'^, 
with E the Ekman number. 



rotation profiles. Some insight into the asymptotic behavior of this mode at small Ekman 
numbers is given by the bottom-right panel in Fig. [9j The damping rate of the mode with 
— 1 is displayed for Ekman numbers varying from 10~^ down to about 2 x 10~^. Runs 



a 



were also performed as a series, the eigenfrequency of the least-damped mode obtained 
at Ekman number E being used as the guess frequency for the run with Ekman number 
E — dE (the spectral resolution is progressively increased with decreasing E). The dashed 

line shows that the damping rate of this particular mode scales proportional to E^'^, a 
scaling that is reminiscent of the solid-body rotation case (Rieutord & Valdettaro 1997). 



The aforementioned E^'"^ scaling of the mode's damping rate is not the only possible 
scaling for shellular rotation in the limit of vanishing viscosities. This is illustrated in the 
top-left panel of Fig. 10 for a particular DT mode having a turning surface that encom- 
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Figure 10. Top-left: evolution of the damping rate (|t|) with decreasing Ekman number for a DT eigenmode 
with shellular rotation, ^(r) oc r°'*, -q = 0.2, which features a detached shear layer tangent to a turning surface 
at the location of the rotation axis. Top-right: contours of kinetic energy for the mode depicted by the rightmost 
open square in the top-left panel [E ~ 8.0 x 1Q~^). The thick white curve shows the attractor of characteristics 
followed by the shear layer, and the dashed curve the location of a turning surface that encompasses the inner 
core. Bottom: cut of the perturbed velocity across the shear layer (at the location depicted by a black line in the 
top-right panel) for the three modes displayed with open square symbols in the top-left panel. 



passes the inner core, and for which the damping rate is found to be roughly independent 
of E. This prop erty is reminiscent of quasi-modes associated with critical layers (see, e.g., 
Le Dizes 2004). However, the physical mechanism leading to such a damping rate must 
be ditterent because the DT mode we consider is axisymmetric, so that no corotation sin- 
gularity exists. Contours of kinetic energy for this new type of quasi-modes are shown in 
the top-right panel of Fig. [lOl where we see that the mode is dominated by a shear layer 



following a wave attractor tangent to the turning surface at the location of the rotation 
axis (note that, since the turning surface encompasses the inner core at r = 0.35/?, this 
mode exists for arbitrarily small inner cores, akin to the mode shown in Fig. |8] for cylin- 
drical calculation). We have carried out additional calculations with a slightly different 
eigenfrequency at given Ekman number, such that the wave attractor no longer hits the 
tangent critical point at the rotation axis. We found that the damping rate is no longer 



independent of E for the same range of Ekman numbers as in the top-left panel of Fig. 10 
It suggests that the tangent critical point plays a prominent role in getting a damping 
rate independent of E. Such scaling demands a very thin detached shear layer, namely 
0(£'^/^)-wide. This is indeed confirmed by the velocity profiles depicted in the bottom 
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Figure 11. Kinetic energy of two axisymmetric D and DT eigenmodes with shellular rotation, for E — 10 . 
Left: D mode with fip « Lll^ref and n{r) oc r"\ Right: DT mode with Q„ ^ 1.90»r cf an d n(r) ex r""'^. Red 



ticks show the values of the inner and outer critical latitudes given by Eqs. (3.281 and (3.291 



panel of Fig. 10 at three Ekman numbers (shown by red squares in the top-left panel of 
Fig. 10 ). Of course, this numerical result does not prove that this vigorously damped mode 
survives at vanishingly small viscosities, and that it is thus a new type of quasi-modes. 
Still, it calls for a detailed analysis of this new type of detached shear layers, which are 
as thin as Ekman boundary layers. 

We now present the results of simulations for two D and DT eigenmodes at smaller 
Ekman number {E = 10~^) and higher spectral resolution (A^^ = 500 or N^ = 580, 
L = 1000). Results are shown in Fig. fiTI a r> ^^^^ „rifi. n ~ i no „ o^^ o ^ ^-i 



11 



A D mode with Up ^ Lllfij-gf and fi oc r 



is shown in the left panel of this figure. A shear lay er is clearly visible that follows an 
attractor of characteristics, calculated through Eq. (3.18) and overplotted by a white 
curve. We note again the excellent agreement betweeii the overall pattern of the viscous 
shear layer, and the paths of characteristics. Because of its finite width, the shear layer 
also excites the critical latitude singularity, on both the inner and out er edg es of the sh ell. 
The location of the inner and outer critical latitudes, given by Eqs. (3.28) and (3.29), is 
depicted by red ticks, and it is clear that the shear layer is tangent to the inner and outer 
edges of the shell at these latitudes. Lastly, the right panel of Fig. ITT] shows a DT mode 



with Up ^ 1.90r2i.ef and f2 oc r~^-^. In this case again, the mode's pattern suggests that 
part of the shear layer approximately follows a wave attractor, shown by a white curve, 
and that part of it is tangent to both the inner and outer edges of the shell at the critical 
latitudes. 

4.2.3. Quasi-periodic orbits of characteristics 

The results presented in the previous paragraphs all feature modes that are associated 
with short-period attractors. These modes are all singular at vanishing viscosity. One may 
however wonder about the existence of regular modes, which exist in the inviscid case, or at 
least modes that display shear layers only at unphysically small viscosities. Such modes are 
necessarily associated with weak attractors or with quasi-periodic orbits of characteristics. 
To assess the occurrence of such orbits of characteristics, we have evaluated the Lyapunov 
exponent of various pairs of characteristics. The Lyapunov exponent measures the rate 
of convergence of two close characteristics after multiple reflections on the inner and 
outer edges of the shell, and on turning surfaces when present. For solid-body rotation, 
the latitudes of two successive reflections on the inner or outer edges can be related 
analytically, and the Lyapunov exponent can therefore be determined semi-analytically 
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(Rieutord et al. 2001). This is no longer the case with differential rotation: Lyapunov 
exponents must be evaluated numerically by integrating paths of characteristics given 
some initial conditions. 

For the symmetry reasons underlined in Sect. |3| paths of characteristics are calculated 
in a meridional quarter-plane with reflections on to the rotation axis (s = 0) and the 
equatorial axis (z = 0). Because turning surfaces may prevent reflections on either axis, we 
deflne the two following Lyapunov exponents associated with reflections on the equatorial 
axis (denoted by Ag) and on the rotation axis (A^): 



A. 



lim 



1 ^' 
-Vln 
V,. z^ 



fc=l 



dsk 



+1 



dsk 



and A, 



lim 



iVz 



— Vln 



k=i 



dzk 



+1 



dzk 



(4.5) 



where ds^ and dz^ are the separation between two characteristics after k reflections on 
the equatorial and rotation axes, respectively. 

A negative Lyapunov exponent means that the two characteristics get progressively 
closer to each other, converging towards an attractor. The more negative this value, the 
faster the convergence (and the stronger the attractor). In contrast, small-amplitude Lya- 
punov exponents may correspond to quasi-periodic orbits of characteristics. Because two 
attractors may co-exist in the shell at a given eigenfrequency, and for averaging purposes, 
we have considered flve pairs of initial conditions and four regions of the shell located 
near the rotation or equatorial axes, and near the shell's surface or the inner core (that 
is, twenty pairs of characteristics). The most negative value amongst the (averaged) Lya- 
punov exponents obtained in the above four regions deflnes what we now refer to the 
Lyapunov exponent, and which we denote by A. If two attractors co-exist in the shell, our 
method will thus retain the one with the shortest period. 

Lyapunov exponents obtained for m = and a few rotation proflles are displayed in 
the top panel of Fig. 12 The range of eigenfrequencies [0.2r2j.ef , 2r2j.ef] is sampled with 
1800 values (f^ref is the angular frequency at the shell's surface). The case with solid-body 
rotation a nd rj = 0.35 is depicted b y a red curve. It can be compared with the semi-analytic 
results of Rieutord et al. (2001), see the left panel of their flgure 6 (where flp/2fij-Qf 
is displayed in x-axis). This comparison shows that our numerical method reproduces 
qualitatively well the ranges of frequencies at which strong attractors and quasi-periodic 
orbits of characteristics exist. Quantitative differences inevitably exist, however, and are 
most signiflcant for attractors with moderate Lyapunov exponents (that is, A in the range 
[—0.5, —0.2]), which our method tends to underestimate. 

Results obtained with shellular rotation and a = 0.5 are overplotted by a black curve in 



the top panel of Fig. 12 For this rotation proflle, inertial modes exist for eigenfrequencies 
|r2p| < 2.2i7i.ef (see, e.g.. Fig. ^ and the transition from D to DT modes occurs at 
\ilp\ ~ l.lSfii-ef (see Eq. 3.25). For Q,p ^ l.lSfiref (D modes only), the shellular case 
features the same attractors as in the case with solid-body rotation. However, attractors 
are shifted toward lower frequencies with shellular rotation as the shell's interior rotates at 
a slower pace in this case. Lyapunov exponents also take very similar values since shellular 
rotation changes the overall shape of these attractors without altering their period. For 
flp ^ 1.18r2j.ef, the shellular case displays many more attractors, which form through 
reflections on to turning surfaces that encompass the inner core. Finally, results obtained 
with the same shellular proflle, but for the full sphere (77 = 0) are also shown by the blue 
curve. The full sphere only supports DT modes in this case. For Qp > lAQj-^^, results are 
identical to those obtained with rj = 0.35, as in both cases the turning surface occupies 
radii r > 0.35i?. Differences for Qp < 1.4r2i-ef are essentially through the strong attractor 
at Qp ~ 0.8f2j.ef, which is related to a turning surface located at r < 0.35-R. 

Although the above results have been obtained with a flxed proflle of shellular rotation, 
the general conclusion that can be drawn is that, just like for solid-body rotation, a 
differentially rotating spherical shell may harbour both short-period attractors and quasi- 
periodic orbits of characteristics. Their occurrence depends on the rotation proflle and 
the sphere's aspect ratio. An example of ?7i = D mode with a weak attractor, obtained 
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Figure 12. Top: Lyapunov exponent obtained for rj — 0.35 and m — 0. Results with solid-body rotation (red 
curve) are compared to shellular rotation with a — 0.5 (black curve). The vertical dashed line displays the 
eigenfrequency that separates D and DT modes for this rotation profile. Results with shellular rotation, a — 0.5 
and for a full sphere (r; = 0) are overplotted in blue. Bottom left: meridional cut of the kinetic energy for an 
m = D mode with a = 0.5 and Qp « 0.92r2rcf, which features a weak attractor (the Lyapunov coefficient is a 
few percent). Bottom-right: meridional cut of the kinetic energy for an rn = 3 D mode with fip ~ — 2.5r2ref and 
a — —0.05. This mode has been obtained assuming asymmetry with respect to the equatorial plane. 



for E = 10^ , a = 0.5 and t] = 0.35, is shown in the bottom-left panel of Fig. 12, and it 
is clear that the kinetic energy for this mode is more homogeneously distributea than in 
the modes illustrated in previous figures. 

While the calculation of Lyapunov coefficients unambiguously shows the existence of 
quasi-periodic orbits of characteristics, it is difficult to use it to prove or disprove the exis- 
tence of reg ular modes for which A should strictly vanish (Dintrans et al. 1999; Rieutord 
). We have looked f or reg ular modes with differential rotation, though. It was 
3y Rieutord et al. ( |2001 ) that toroidal modes, with Doppler-shifted eigenfre- 

quencies of the form Qp = 2Q/{m + 1) {m is the azimuthal wavenumber) are the only 
regular modes in a rigidly-rotating spherical shell. We have followed the toroidal mode 



et al.\\2001 
suggested 
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Figure 13. Top: meridional cut olthe viscous dissipation (left quadrant) and kinetic energy (right quadrant) for a 
D mode with rra = 2, eigenfrequency Q,p ~ —2.3Q,Tei, E = 2x 10~® and cylindrical rotation profile Q,{s) — 1 + 0. 4s^. 
A corotation resonance (critical cylinder) is located at s ~ 0.6. Bottom: spectral content of the radial (it) and 
orthoradial (w) components of the velocity field for this mode. Chebyshev and spherical harmonic coefficients are 
shown in the left and right panels, respectively. 



with TTi = 3 by slowly increasing the differential rotation rate (for shellular rotation). As 
shown in the bottom-right panel of Fig. 12, even for values of a as small as a few percent, 
the mode's regularity fades away. Although most of the kinetic energy remains contained 
in a smooth structure reminiscent to that of the toroidal mode, ray structures develop, in 
particular from the inner critical latitude. Therefore, inertial oscillations with differential 
rotation seem to be all singular in the lim it of vanishing Ekman numbers. A simila r result 



was obtained for gravito-inertial modes (Friedlander 1982 Dintrans et al. 1999), which 



also satisfy a PDE of mixed type, analogous to Eq. (tTir 



4.3. Non-axisymmetric eigenmodes: corotation resonances 
We briefly examine in this section the behavior of viscous, non-axisymmetric inertial waves 
when their frequency in the rotating frame (f2p) vanishes inside the shell. Such corotation 
resonances will be referred to as critical cylinders for cylindrical rotation, and as critical 
layers for shellular rotation. We present results of linear numerical calculations, which 
we compare to the an alysis of Sectio n |3} Cylindrical and shellular rotation profiles are 
considered in Sections |4.3.1| and |4.3.2[ respectively. 



4.3.1. Cylindrical rotation 

An example of m = 2 D eigenmode with a critical cylinder is displayed in Fig. [13 
This mode is obtained for e 

Nr = 880 X L 



0.4, E' = 2 X 10 ^, and the very high spectral resolution 
2400. Its eigenfrequency Qp ~ — 2.3r2i-ef (a negative Qp is required for 



Inertial waves in a differentially rotating spherical shell 



25 




O.D 0.2 0.4 0.6 O.S 1.0 

Nr=600 L=1400 M=2 E=5,0x10-'ti=0.350e=0.4O0 




200 400 600 800 1000 1200 1400 
I 



Figure 14. Left: meridional cut of the kinetic energy for a D mode with m — 2, eigenfrequency fip ~ — 2.20r2rcf , 
E = 5 X 10~^, and e = 0.4. A corotation resonance is located at s ~ 0.5. A few paths of characteristics are 
overplotted by white curves. Right: spectral content of the radial (u) and orthoradial (w) components of the 
velocity field for this mode. Chebyshev and spherical harmonic coefficients are shown in the top and bottom 
sub-panels, respectively. 

an m > eigenmode to have a corotation resonance within the shell). The top panel of 
Fig. 13 shows a meridional cut of the viscous dissipation (left quadrant) and of the kinetic 
energy (right quadrant) for this mode. The mode's spectral content is also displayed in the 
bottom panel of Fig. [13] to highli ght t he good convergence in spectral resolution obtained 
for this mode. AccorcWg to Eq. (3.13), a critical cylinder is located at Scrit ~ 0.61-R. Near 
this location, the mode's kinetic energy and viscous dissipation are discontinuous. Upon 
approaching corotation from s > Sent, the mode's pattern in a meridional qu arter-p lane 
becomes more and more parallel to the rotation axis, and according to Eq. (3.16), the 
component of the group velocity perpendicular to the critical cylinder should eventually 
cancel out at corotation, formally preventing the mode from crossing the critical cylin- 
der. Our simulation indicates that, in a moderately viscous fluid, waves get dramatically 
absorbed upon crossing corotation, but keep on propagating at s < Scrit with a much 
reduced amplitude. 

We point out that the fraction of the mode's amphtude that is effectively transmitted 
across corotation depends on the mode of interest. A similar m = 2 mode with eigen- 
frequency and Ekman number slightly increased to fln^ — 2.20firef and E = 5 x 10~ , 
respectively, shows no significant transmission in Fig. [Ml as if the incoming wave was fully 



absorbed at corotation. The contour plot in the left panel of this figure gives an upper 
limit for the kinetic energy transmission factor ~ 10~^. Note that in the example shown 
in Fig. 14, linear inertial waves propagate from the inside outwards, whereas Fig. 13 shows 
an opposite direction of propagation (from the outside inwards). It is likely due to a mode 
selection, with different modes existing on each side of the corotation radius. 

4.3.2. Shellular rotation 

An example oi m = —6 eigenmode with a critical layer is depicted in Fig. 15 For this 
mode, E = Q X 10~^, VL{r) oc r~^'^, and the spectral resolution is Nr = 500 x L = 1401. 
This mode's eigenfrequency Qp ~ T-llfij-ef (this time, a positive Qp is required for an 
m < eigenmode to have a corotation resonance within the shell). Contours of the 
kinetic energy in a meridional quarter-plane are shown together with the mode's spectral 
content. The critical layer is depicted by a white solid curve in the left panel, and we see 
that, in contrast to cylindrical rotation, waves cross the critical layer without experiencing 
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Figure 15. Left: meridional cut of the kinetic energy for a growing eigenmode with iJ = 6xlO '■',cr = — 0.2, 
7.1inrcf and m = —6. The location of the critic al lay er is depicted by a solid curve, along which a white tick 



i Lp 

shows the local slope of characteristics given by Eq. (3.341. Red ticks show the inner and outer critical latitudes 
Right: spectral content of the radial (the quantity v — xuis displayed here) and orthoradial (w) components of the 
velocity field for this mode. Chebyshev and spherical harmonic coefficients are in the top and bottom sub-panels, 
respectively. 



significant absorption. The local analysis carried out for the inviscid case in Section |3.2.3 
shows that (linear) inertial waves may indeed cross a critical layer with a (formally) 
divergent group velocity at the critical layer. This occurs ii\kz\ — )■ (paths of characteristic 
then being locally parallel to the rotation axis), or if i3 — > (paths of characteristics having 
a finite slope at the critical layer). The white tick displayed on top of the critical layer in 
the l eft pa nel of Fig. 15 shows the local slope of characteristics if i3 = 0, which is given by 
Eq. (3.34). The large slope obtained in this case (due to the small value of \a\) prevents a 
clear distinction with the case l^^l — )• 0. The local analysis also shows the possibility that 
the group velocity formally vanishes at a critical layer (when \ks\ — > oo locally), much 
like in the case of cylindrical rotation. However, we have not found such situation in our 
numerical exploration of modes featuring a critical layer. 

The most salient feature of the mode shown in Fig. [15] is that it is found to be unstable: 
its growth rate, r, is positive. We have examined in some more details this particular mode 
by varying the Ekman number, keeping the same spectral resolution. We have varied E in 
the range [4 x 10~^ — 3.2 x 10~®]. Results are shown in Fig. 16, where the mode's growth 
rate is depicted as a function of the Ekman number. There is a transition occurring near 
Ec ~ 10~®: r < (the mode is decaying) for E > Ec, whereas r > (the mode is 
growing) for E < Ec- The mode's growth rate is actually well approximated by a linear 
function of E near this transition: r = a{E — Ec), with a < 0. Although not displayed 
here, the kinetic energy obtained in the simulations with E > 10~® is very similar to that 
of the left panel in Fig. 15 , except for the thick shear layer on each side of corotation, 
which is found to progressively disappear with increasing E. The very good convergence 
in spectral resolution obtained in these calculations (as is highlighted by the right panels 
in Fig. 15), and the fact we have been able to successfully follow this particular mode with 
decreasing Ekman number (see Fig. 16 ) makes us confident that the positive growth rate 



that we obtain at small Ekman nunTbers is physical, and not a numerical artifact arising 
from truncature or round-off errors. It is temp ting to interpret this result based on the 
local dispersion relation obtained at Eq. ( 3.30[ ). It is possible indeed that the instability 
we have obtained in our simulations is related to the Goldreich-Schubert-Fricke instability 
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Figure 16. Growth rate for the same eigenmode as in Fig. 15 {m = —6, a — —0.2, Q,p ~ 7.11Jlrcf) with varying 
the Ekman number E. The mode becomes unstable for E < 10~*. 

in the limit of an incompressible inviscid fluid. The physical mechanism that drives the 
instability in our model requires a speciflc investigation that is postponed for future work. 

5. Astrophysical implications and concluding remarks 

The late evolution of the orbital elements and spin of short-period extrasolar planets, 
or of close binary-star systems, is primarily driven by tidal interactions. The tidal torque 
acting on a gravitationally perturbed body encompasses the equilibrium torque, which 
arises from the quasi-hydrostatic tidal bulge setting around the body, and the dynamical 
torque, due to tidally forced internal waves propagating inside the body. Inertial waves 
may propagate in convective regions of stars and planets, where the Coriolis force is 
the main restoring force. Gravito-inertial waves may propagate in radiative zones, where 
buoyancy and rotation act as restoring forces. The magnitude of the dynamical torque is 
directly related to how much angular momentum is deposited by internal waves, which 
has a complex dependency upon the forcing potential (amplitude and frequency), and the 



body's internal p roperties (see, e.g., Ogilvie & Lin 2004, 2007 Ogilvie 2009 Rieutord &; 



Valdettaro 2010). In many dynamical studies of short-period exoplanets, the tidal torque 



ot the central star is parametrized by a single, dimensionless tidal quality factor Q, whi ch 



should be thought of, and determined as, a frequency- averaged quantity (Ogilvie 2009). 

Prior to an investigation of tidally-forced inertial waves in differentially rotating stars 
and planets, we have examined the propagation properties of linear inertial waves in 
a rotating fluid body contained in a spherical shell, taking into account the effects of 
differential rotation. We have used a simple model that consists of an incompressible 
homogeneous fluid, with dissipation occurring through viscous forces. This fluid model is 
an intentionally simplifled model for the convective zones of stars and planets. We have 
considered two particular rotation proflles: a cylindrical rotation proflle, and a shellular 
rotation proflle. For the sake of simplicity, we have assumed the background flow to be 
in equilibrium with differential rotation. In the case of shellular rotation, this equilibrium 
would require additional forces that we have not included. In this respect, the relaxation 
of the incompressible assumption, or the generalization to stratifled fluid bodies, will 
certainly be an improvement on our model. 

Upon linearization of the governing fluid equations, we have revisited the propagation 
properties of paths of characteristics with differential rotation, in the inviscid limit. Paths 
of characteristics are shown to obey a second-order PDE of mixed type. Two kinds of 
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inertial modes are thus found to exist with differential rotation: D modes, for which the 
hyperbohc domain covers the entire shell, and DT modes, for which the hyperbolic do- 
main covers part of the shell. We have analyzed the occurrence of axisymmetric {m = 0) 
and non-axisymmetric D and DT eigenmodes for various differential rotation rates and 
eigenfrequencies. The presence of turning surfaces for DT modes extends the range of fre- 
quencies at which inertial modes of oscillation may exist compared to solid-body rotation. 
We have also discussed the presence of critical latitudes and corotation resonances with 
differential rotation. 

The analytic results based on the dynamics of characteristics are then compared to 
the results of numerical simulations of linear inertial waves propagating in a viscous fluid 
contained in a spherical shell. High-resolution calculations, based on a spectral method, 
have been performed to properly describe the structure of the thin shear layers arising in 
the physically interesting regime of small Ekman numbers (small viscosities). The Ekman 
numbers used in our simulations typically range from a few 10~^ to 10~^, values that are, 
however, much larger than what may be expected in stars and planets (Ekman numbers 
between 10~^^ and 10~^^). The modes that we have primarily focused on are singular 
modes with a shear layer converging towards a short-period wave attractor. As in the solid- 
body rotation case, the pattern of these shear layers is well reproduced by the dynamics of 
characteristics under the short- wavelength approximation. As an interesting side result, 
we have shown that reflexions on a turning surface, and on the body's outer radius, 
may form singular modes with short-period wave attractors in the whole sphere (that is, 
in the absence of an inner core). The damping rate of singular modes with differential 
rotation is also assessed in the limit of small Ekman numbers. Our results of calculations 
indicate that an E ' scaling is to be expected for most D and DT eigenmodes, much 
similar to the case of solid-body rotation. Still, the case of a DT mode with damping 
rate approximately independent of the Ekman number is uncovered. In addition, we have 
carried out a qualitative analysis of the occurrence of modes with quasi-periodic orbits of 
characteristics for shellular rotation. 

We have also briefly examined the behavior of linear inertial waves with a corotation 
resonance within the shell (where the wave's frequency vanishes in the frame rotating with 
the fluid). For cylindrical rotation, our simulations show that inertial waves may cross a 
corotation resonance, but experience a very large absorption at this lo cation, much like 



inertial waves propagating in differentially rotating astrophysical discs (Latter & Balbus 
12009). For shellular rotation, however, we have found that waves may propagate across a 
corotation resonance (critical layer) without visible absorption. Such modes are found to 
be unstable for small enough viscosities. Clearly, these findings appeal for a confirmation 
with an improved fluid model, where shellular rotation is balanced for by appropriate 
body-forces, and where non-linear effects are taken into account. 

Many simplifying assumptions in the fluid model we have considered deserve further 
investigation to make it more applicable to the fluid flows inside stars or planets. The 
incompressible and homogeneous assumptions are obvious examples, and naturally a full 
model including appropriate modeling for the convective and radiative zones would be de- 
sirable. Since inertial waves propagate in convective regions of stars and planets, it would 
be interesting to study the impact of convective turbulence on the propagation and dis- 
sipation properties of inertial waves (the formation of large-scale attractors for instance). 
The nature and quality of wave reflections at the shell's boundaries should be analyzed, 
as it is still unclear for example how inertial waves reflect on the upper atmosphere of 
stars and planets. 
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Appendix A. 

This appendix derives the partial differential equation satisfied by the perturbed pres- 
sure, p, in the inviscid limit. It also derives the wave dispersion relation, the phase and 
the group vel ociti es under the short-wavelength approximation. In cylindrical coordinates 
{s,ip,z), Eq. (2.1) becomes 



ifipU + 2i7ez X u + s(u ■ Vi7)e^ = Vp, 



(Al) 



with Q the background angular frequency and Qp the Doppler-shifted frequency. The 
cylindrical components of the linearized velocity vector, {us,Uip,Uz), can be expressed as 
a function of p and of its partial derivatives: 
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and where 



A, = ^^(s2fi) and Az = —^{s^n). 



(A 5) 
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For solid-body rotation or cylindrical rotation, A^ = and As = 4^2^. Note the special 
case where Qp = for which a corotation resonance exists in the shell, and that where 
As = fip, which is analogous to Lindblad resonances in compressible fluids with cylindrical 
rotation such as astrophysical discs. 



Taking the divergence of Eq. (Al) and using Eqs. (A 2) to (A 4), we find, after some 
algebra, that p satisfies: 
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which reduces to Eq. (3.1 ) when only second-order terms are retained (the short-wavelength 
approximation) . 

We now examine the propagation prop erties of local waves in a meridional plane. As- 
suming p oc exp{i{ksS + k^z}), Eq. (A 6) yields the following wave dispersion relation in 
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the short-wavelength approximation: 



nt 
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(A 7) 



where ||k|| = ykl + k^. For sohd-body rotation, Eq. (A 7) reduces to the well-known 
dispersion relation Ctp = ±2Qkz 



[see, e.g. 



Greenspan 1969). In the bracket term of 



Eq. (A 7), while Ag > i s req uired to satisfy Rayleigh's criterion, the sign of Azks/kz is 
not known a priori. Eq. (A 7) therefore indicates that an instability may occur if Ag — 
Azkg/kz < 0. It corresponds in fact to the Goldreich-Schubert-Frick e instability for an 
incompressible inviscid fluid (Goldreich & Schubert 1967 Fricke 1968), in whic h case the 
inequality As — A^ k.^/kz < is identical to that in equation 33T niisprinted) of 
& Schubert (1967) (see also Nelson, Gressel & Umurhan (2012)). Assuming t 
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stable against the Goldreich-iSchubert-l^'ricke instability, we introduce 
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and Eq. (A 7) translates into ^p = ±i3|A;2|/||k||. In the frame rotating with the fluid, the 
phase velocity, Vp = i7pk/||kp, satisfies 



while the group velocity, Vg = dflp/dk, is given by 
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Note that k-Vg = 0: the group velocity is perpendicular to the phase velocity for differential 
rotation as for solid-body rotation. 
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